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A novel approach that enables the study of parallel transport in magnetized plasmas is presented. 
The method applies to general magnetic fields with local or nonlocal parallel closures. Temperature 
flattening in magnetic islands is accurately computed. For a wave number k, the fattening time 
scales as xw''' ~ k~" where x is the parallel diffusivity, and a = 1 (a = 2) for non-local (local) 
transport. The fractal structure of the devil staircase temperature radial proflle in weakly chaotic 
fields is resolved. In fully chaotic fields, the temperature exhibits self-similar evolution of the form 

T = {x\\t)~^'^L (X||*) ' '^'Z' 1 where Sip is a radial coordinate. In the local case, / is Gaussian 

and the scaling is sub-diffusive, 7 — 1/2. In the non-local case, / decays algebraically, L{ri) ~ 77"^, 
and the scaling is diffusive, 7=1. 
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The study of transport in magnetized plasmas is a 
problem of fundamental interest in controlled fusion, 
space plasmas, and astrophysics research. Three issues 
make this problem particularly challenging: (i) The ex- 
treme anisotropy between the parallel (i.e., along the 
magnetic field), X|h ^^^ ^^^ perpendicular, X-L, conduc- 
tivities (x||/X-L i^ay exceed 10^° in fusion plasmas); (ii) 
Magnetic field lines chaos which in general complicates 
(and may preclude) the construction of magnetic field 
line coordinates; and (iii) Nonlocal parallel transport in 
the limit of small coUisionality. As a result of these chal- 
lenges, standard finite difference and finite elements nu- 
merical methods staffer from a number of ailments. Chief 
among them are the pollution of perpendicular dynamics 
due to truncation errors in the discrete representation of 
the parallel heat flux, and the lack of a discrete maximum 
principle (to enforce temperature positivity). 

Despite the severity of these issues, recent studies have 
succeeded in partially addressing some of them, and im- 
portant progress has been made in the study of parallel 
transport. Reference HE] discussed a finite element numer- 
ical implementation of nonlocal heat transport with ap- 
plications including temperature flattening across mag- 
netic islands and tokamak disruptions. The use of high- 
order discretizations has been shown to mitigate numer- 
ical pollution of the perpendicular dynamics in finite- 
difference [HIS] and finite-element methods [31I3]. A max- 
imum principle has been shown to be enforceable by the 
use of limiters at the discrete level in finite differences [5] , 
and finite elements [6], albeit with the effect of render- 
ing both spatial discretizations formally first-order accu- 
rate. In Ref. [7j a second order finite-difference iterative 
Krylov method was used to find the steady state solu- 
tion of the heat transport equation in a weakly chaotic 
magnetic field. 

Motivated by the strong anisotropy typically encoun- 
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tered in magnetized plasmas (xii/X-L ^ 10^*^) i "wo study 
parallel heat transport in the extreme anisotropic regime 
(X||/X-L ~^ 00) ■ Ii^ addition to the previously mentioned 
numerical difficulties, this regime presents potentially in- 
surmountable issues regarding the algorithmic inversion 
of the discretized transport equation. In particular, the 
potentially degenerate null space of the parallel transport 
operator might render the use of state-of-the-art scalable 
iterative inversion methods (e.g. multigrid methods) im- 
practical. To overcome these numerical challenges, we 
present a novel Lagrangian Green's function approach. 
The proposed method bypasses the need to discretize 
and invert the transport operators on a grid and allows 
the integration of the parallel transport equation without 
perpendicular pollution, while preserving the positivity 
of the temperature field. The method is applicable to lo- 
cal and non-local transport in integrable, weakly chaotic, 
and fully chaotic magnetic fields. 

Beyond the integration method, this Letter presents 
novel physics results. The difference between local and 
non-local parallel temperature mixing and flattening in- 
side magnetic islands is studied. Of particular interest 
is the dependence of the temperature relaxation time 
on the wave number of the temperature perturbation. 
Understanding this problem is of significant interest to 
assess the impact of magnetic islands on magnetic con- 
finement. Parallel transport is also studied in the case 
of weakly chaotic and fully chaotic 3-D magnetic fields. 
In the case of weakly chaotic fields, the fractal Hamil- 
tonian structure of the magnetic islands implies that, as 
the ratio X\\/x± increases, the radial temperature pro- 
file approaches a devil-staircase. Going beyond previous 
studies [7], we unveil the fractal structure of the devil- 
staircase in the previously unaccessible X-L = regime. 
These results open the possibility of a deeper understand- 
ing of the role of cantori which have been observed to act 
as partial transport barriers in numerical studies VT\ and 
experiments [& . Another problem of considerable inter- 
est in fusion and astrophysical plasmas is the understand- 
ing of electron heat transport in fully chaotic fields. The 



case of local parallel transport has been extensively stud- 
ied since the pioneering work in Refs. [9l[l0]. However, 
most studies have restricted their attention to local trans- 
port. Here, in addition of the study of the local transport, 
we present novel results on the self-similar, non-Gaussian 
spatio-temporal scaling of radial heat transport in the 
non-local case. 

Our starting point is the heat transport equation in a 
constant density plasma 



dtT 



Vq, 



(1) 



where q is the heat flux. For local transport, in the limit 
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(2) 



where b = B/|S|. Substituting ([2|) into ([ij), we get 

dtT^~dsq\\. qi\=-x\\dsT, (3) 

where ds = b • V is the derivative along the field line, 
and we have assumed the tokamak ordering dghiB « 0. 
In the case of non-local transport, following Ref. |TT], we 
consider 



1\\ 



Ax 



T{s + z)-T{s- z) 



dz . 



(4) 



with 1 < a < 2. The case a = 1 reduces to the col- 
lisionless heat flux characterized by the free-streaming 
of electrons along magnetic field lines HH IS] , and the 
case a — 2 reduces to the local diffusive case. The 
model in Eq. Q allows the interpolation between these 
two regimes. In particular, Eq. Q allows the incorpora- 
tion of parallel transport processes with underlying non- 
Gaussian (Levy a-stable) stochastic processes |11| . 

Substituting Eq. Q into Eq. ([T|, and assuming X|| con- 
stant, the transport equation canbe written in a compact 
form as 



dtT = x\\d^,lT 



(5) 



where the operator cJP, denotes the symmetric fractional 
derivative of order a, defined as J^ df^^T = —\k\"T. 

Here, J-" denotes the Fourier transform [TT]. As expected, 
in the limit a = 2, Eq. ^ reduces to the diffusion equa- 
tion in Eq. ^. 

The proposed method is based on the Green's function 
solution of Eq. ([5| along the Lagrangian trajectories of 
the magnetic field. The unique magnetic field line trajec- 
tory (parametrized by the arc length) that goes through 
a point Tp is given by the solution of 



-=b 

ds 



r(s = 0) 



Lp. 



(6) 



Thus, given an initial condition in the whole domain, 
r(r,i = 0), the initial condition along the field line is 





FIG. 1: Temperature mixing and homogenization in a mag- 
netic island. Panel (a) shows the initial condition To = 
2-|-cos(4^), and panel (6) the final relaxed state, Tf = F{Az), 
where A^ is the magnetic potential, for both local and non- 
local transport. Panels (c) and (d) show the solution at 
X||i = 0.6 in the local and non-local cases respectively. 



Tq{s) — T{r{s),t — 0) and the temperature at r = Tp at 
time t is given by 



T{rp,t)= I \oHs)]Gc.{s,t)ds, 

J si 



(7) 



where Gq is the Green's function of Eq. (Is]). For un- 
bounded field lines, (si,S2) — (— cxi,cxd) and the Green's 
function is given by 



G„(s,i) 



I 

2^ 



-Xt|fc|°-ifcs^^ 



For a = 2, Eq. (|8|) gives the Gaussian distribution 



G2(s,i) = ^(x|ii) 



-1/2 



exp 



s 

4^ 



(8) 



(9) 



and in the non-local free streaming case, a = I, it gives 
the Cauchy distribution 



Gi(s,t) 



(Xil^)"^ 



^ + {s/x\\ty 



(10) 



For general a, Gq ~ (x||t)^"^^"ia,o[(x||i) "^^"■s] where 
Lq.o is the symmetric a-stable Levy distribution, see for 
example Ref. [14. In the case of closed (periodic) field 
lines, the integration domain in Eq. ([7|) is a finite inter- 
val (si,S2) with r(si) = r(s2), and one has to use the 
periodic Green's function, G^, obtained by mapping the 
unbounded Green's function, G^, into the periodic do- 
main. 




FIG. 2: Poincare plot of the weakly chaotic magnetic field 
used in solution of the parallel heat transport equation shown 
in Fig. [I 
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FIG. 3: Radial temperature profile of the time-asymptotic 
solution of the parallel heat transport equation for the weakly 
chaotic magnetic field in Fig. [2] The zooms in the successive 
panels unveil the fractal structure of the devil-staircase profile. 



At this point it is important to indicate a fundamen- 
tal difference between the work in Refs [U [T3] and the 
method proposed here. In Refs. [IIII3], the flux is calcu- 
lated throughout the computational domain by integrat- 
ing along the field lines. However, in these references, 
once the flux is computed, the Lagrangian approach is 
abandoned and the flux is mapped to Gaussian quadra- 
ture points for the flnite-element standard integration of 
the temperature evolution equation on a grid. On the 
other hand, the method proposed here is fully Lagrangian 
and completely bypasses the use of finite differences or fi- 
nite elements integration schemes, circumventing the nu- 
merical limitations discussed before. Another clear ad- 
vantage of the use of the Green's function is that the 



solution at an arbitrary time i = t/ is obtained directly 
from the integral in Eq. ^\. This is to be contrasted 
with the standard explicit and implicit time discretiza- 
tion schemes that require the knowledge of the solution 
at previous times t < tf before the solution at t — tf can 
be constructed. A similar decoupling occurs in space. 
According to Eq. (l7|, the evaluation of T at r^ is fully 
decoupled from the solution at r 7^ rp, yielding naturally 
to a massively parallel approach. 

In what follows, we use the Lagrangian Green's func- 
tion method to study local and non-local parallel trans- 
port in magnetic fields of increasing levels of complexity. 
We start with the study of transport in magnetic islands, 
which we model with a 2-D magnetic field with vector po- 
tential Az = if/)^ + ecos6', where ip and 9 are radial and 
angular coordinates respectively. Of particular interest 
is the difference in the temperature relaxation properties 
in the presence of local and nonlocal transport. Figure 1 
shows the solution of the parallel transport of an initial 
condition of the form Tq = 2 + cos(fc6') with fc = 4. As ex- 
pected, for large enough times, the temperature relaxes 
to a unique solution of the form Tf = F{Az) in both the 
local and non-local cases. However, the dynamic process 
leading to the final relaxed state is different. As the fig- 
ure shows, temperature mixing in the non-local case is 
less homogeneous and slower than in the local case. In 
general, from Eq. (l8| it follows that the homogenization 
time of a temperature perturbation with mode number 
k scales as xw't ~ k~°' with a — 2 {a = 1) in the case of 
local (non-local) transport. 

For the study of transport in 3-D magnetic fields, we 
assume cylindrical geometry periodic in z with period 
L = 2nR where i? is a constant. The magnetic field 
consists of a perturbed screw-pinch of the form 

B = {rB/\)/[l + {rl\f]eg + B^e, + Bi(r, 9, z) , (11) 

with B, A, and Bq constants. The magnetic potential of 
the perturbation, Bi = V x A^e^, consists of a superpo- 
sition of modes 



A.{r,9,z) 



Z^ 



Amn{r) cos {m9 - nz/R + Cmn) , (12) 



with 



^n 



ea{r) 



exp 
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V2a 



ro 
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(13) 

For each {m, n), the values of r*and tq are chosen so that 
the safety factor satisfies q{r.^,) = m/n and dA„in/dr(r = 
r*) = 0. The prefactor (r/r^,)"'^ is included to guarantee 
the regularity of the radial eigenfunction near the origin, 
r ~ 0. The function, a{r) = [1 - tanh[(r - l)/l]]/2, is 
introduced to guarantee the vanishing of the perturbation 
for r ^ 1 and the existence of well-defined flux surfaces 
at the plasma boundary. 

In the study of transport in weakly chaotic flelds, only 
two modes were included. As the Poincare plot in Fig- 
ure [2] shows, the magnetic fleld in this case exhibits a 
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FIG. 4: Self-similar spatio-temporal evolution of the radial 
temperature profile for local, a = 2, transport in a fully 
chaotic magnetic field. In this case the scaling function is 
Gaussian and the scaling exponent is 7 = 1/2. 
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FIG. 5: Self-similar spatio-temporal evolution of the radial 
temperature profile for non-local, a = 1, transport in a fully 
chaotic magnetic field. In this case the scaling function is 
strongly-non-Gaussian and, as the dashed-line fit shows, ex- 
hibits algebraic decay of the form L ~ r;"'^. The scaling ex- 
ponent is 7 = 1. 



rich fractal-like structure resulting from the existence of 
higher-order resonances. Figure [3] shows the time asymp- 
totic, final radial temperature solution along the = 2.7 
horizontal line in Fig.j2J corresponding to the initial con- 
dition To = 1 — 2i?^. As expected, the same state 
is reached in the presence of local and non-local trans- 
port. However, as mentioned before, the relaxation time 
is longer in the non-local case. The fractal structure of 
the magnetic field gives rise to a devil-staircase temper- 
ature profile in which horizontal sections (resulting from 
high-order resonances) are mixed with vertical sections 
(resulting from KAM invariant circles and Cantori). 

To study transport in a fully chaotic magnetic field we 
considered a set of 21 strongly overlapping modes. In this 
case, the Poincare plot (not shown) is fully hyperbolic 



and does not exhibit any structure. The initial condition 
consists of a narrow "cylindrical shell" of the form Tq — 
exp[-i?2(V; - '(/'o)/CT]^ with ^0 = 0.18 and a = 0.05. 
Figures |4] and [5] show the time evolution of the radial 
profile of the temperature averaged in 6 and z in the local 
and the non-local free streaming cases. In both cases, the 
temperature exhibits an asymptotic self-similar evolution 
of the form 



{T)oA^,t) = ixtr''^^L{rj)., 



(14) 



where the similarity variable is defined as rj = {ip— < 
-0 >)/{xtV^'^ with 7 the scaling exponent. From here, it 
follows that the second moment scales as < ijP' >^ V . 
As Fig. |4| shows, in the local transport case the scaling 
function is approximately Gaussian, L ^ e~^ ''^. Consis- 
tent with Ref. [10], in this case 7 w 1/2. This subdiffu- 
sive scaling results from the combination of the diffusive 
transport along the field line and the quasilincar diffusion 
of the field line itself due to the magnetic field line chaos. 
On the other hand, in the nonlocal case, the scaling func- 
tion is strongly non-Gaussian and exhibits an algebraic 
decay of the form L ~ 77"^. Interestingly, in the non- 
local case 7^1. This diffusive scaling results from the 
coupling of the quasilinear diffusion of the field line due 
to chaos and the free-streaming transport along the field 
line. 

Summarizing, this Letter presented a Lagrangian 
Green's function method for the accurate and efficient 
computation of purely parallel (xj_ = 0) local and non- 
local transport in arbitrary magnetic fields with constant 
X|| • Because of the parallel nature of the Lagrangian cal- 
culation, the formulation naturally leads to a massively 
parallel implementation, suitable for today's supercom- 
puters. The method was applied to compute temperature 
mixing in magnetic islands, for which the ratio of the 
non-local free-streaming, {x\\T)ni, and diffusive, (x||T)rf, 
transport relaxation times scales as {x\\''')ni/{x\\T)d ~ k. 
Radial transport in 3-D magnetic fields in cylindrical ge- 
ometry was also studied. In the case of weakly chaotic 
fields, the fractal structure of the devil staircase of the ra- 
dial temperature profile was resolved. In the case of fully 
chaotic fields, the radial transport exhibited self-similar 
spatio-temporal behavior. Contrary to the well-known 
diffusive case, the non-local case exhibits a non-Gaussian, 
algebraic decaying temperature profile with scaling expo- 
nent 7 = 1. Current work includes the implementation of 
the Lagrangian Green's function method as integral part 
of a more general multiscale framework that incorporates 
finite x_L, non-constant x\\i and heat sources. 
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